use rand::Rng;


pub fn addmod(a: u128, b: u128, p: u128) -> u128 {
    let (a, c) = a.overflowing_add(b);
    if c {
        (a + (!p) + 1) % p
    } else {
        a % p
    }
}

pub fn mulmod(mut a: u128, mut b: u128, p: u128) -> u128 {
    let mut res = 0;
    while b > 0 {
        if (b & 1) == 1 {
            res = addmod(res, a, p);
        }
        a = addmod(a, a, p);
        b >>= 1;
    }
    res % p
}

fn power(mut a: u128, mut b: u128, p: u128) -> u128 {
    let mut res = 1;
    while b > 0 {
        if (b & 1) == 1 {
            res = mulmod(res, a, p);
        }
        (a, b) = (mulmod(a, a, p), b >> 1);
    }
    res % p
}

fn break2(mut x: u128) -> (u128, u128) {
    let mut cnt = 0;
    while (x & 1) == 0 {
        x >>= 1;
        cnt += 1;
    }
    (x, cnt)
}

fn mr(n: u128) -> bool {
    if n < 2 || n % 2 == 0 {
        return n == 2;
    }
    let (d, c) = break2(n - 1);
    let arr =vec![2, 3, 7, 61, 24251, 4840261];
    for a in arr {
        if a >= n {
            continue;
        }
        if n % a == 0 {
            return false;
        }
        let mut t = power(a, d, n);
        if t == 1 || t == n-1 {
            continue; // pass
        }
        for _ in 1..c {
            t = mulmod(t, t, n);
            if t == 1 { // fail
                return false;
            }
            if t == n - 1 { // pass
                break;
            }
        }
        if t != n - 1 { // fail
            return false;
        }
    }
    true
}

fn gcd(x: u128, y: u128) -> u128 {
    if y == 0 {
        x
    } else {
        gcd(y, x % y)
    }
}

fn dist(x: u128, y: u128) -> u128 {
    if x >= y {
        x - y
    } else {
        y - x
    }
}

fn pollards(n: u128) -> u128 {
    let mut rng = rand::prelude::thread_rng();
    let mut x = 0;
    let mut y = 0;
    let c = rng.gen_range(1..=n);
    let mut k = 1;
    loop {
        let mut t = 1;
        let mut q = 1;
        for i in 0..k {
            x = addmod(mulmod(x, x, n), c, n);
            q = mulmod(q, dist(x, y), n);
            if (i & 0x7F) == 0 {
                t = gcd(q, n);
                if t > 1 {
                    break;
                }
            }
        }
        if t > 1 {
            break t;
        }
        t = gcd(q, n);
        if t > 1 {
            break t;
        }
        k *= 2;
        y = x;
    }
}

pub fn find_max_prime_factor(mut number: u128) -> u128 {
    if mr(number) || number == 1 {
        number
    } else {
        let mut tmp = number;
        while tmp == number {
            tmp = pollards(number);
        }
        while number % tmp == 0 {
            number /= tmp;
        }
        let ans1 = find_max_prime_factor(number);
        let ans2 = find_max_prime_factor(tmp);
        std::cmp::max(ans1, ans2)
    }
}
